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GEOMETRICALLY  NONLINEAR  TRANSIENT  ANALYSIS  OF 
LAMINATED  COMPOSITE  PLATES 


J.  N.  Reddy t 

Virginia  Polytechnic  Institute  and  State  University 
Blacksburg,  Virginia  24061 

Abstract 

Forced  motions  of  laminated  composite  plates  are  Investigated  using 
a  finite  element  that  accounts  for  the  transverse  shear  strains,  rotary 
inertia,  and  large  rotations  (In  the  von  Karman  sense),  the  present 
results  when  specialized  for  isotropic  plates  are  found  to  be  In  good 
agreement  with  those  available  In  the  literature.  Numerical  results  of 
the  nonlinear  analysis  of  composite  plates  are  presented  showing  the 
effects  of  plate  thickness,  lamination  scheme,  boundary  conditions,  and 
loading  on  the  deflections  and  stresses.  The  new  results  for  composite 
plates  should  serve  as  bench  marks  for  future  Investigations. 
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extenslonal ,  flexural -extenslonal  and  flexural 
stiffnesses  (1,j  =  1,2,6). 
plate  planform  dimensions  of  x,y  directions, 
parameters  in  the  time  approximation  (see  eqn.(15)) 
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rotatory  Inertia  coefficient  per  unit  midplane  area  of 
layer 

shear  correction  coefficients  associated  with  the  yz 
and  xz  planes,  respectively  (1  =  4,5) 
stress  couple  and  stress  resultant,  respectively 
(1=1. 2, 6) 

laminate  normal  Inertia  coefficient  per  unit  midplane 
area 

shear  stress  resultant  (1  =  1,2) 
plane-stress  reduced  stiffness  coefficients  (i,j  = 
1,2,6) 

laminate  rotatory-normal  coupling  inertia  coefficient 
per  unit  midplane  area 

displacement  components  In  x,y,z  directions, 
respectively 

nodal  values  of  displacements  u,v,w  (i  =  l,2,...,n) 
position  coordinates  in  cartesian  system 
column  vector  of  generalized  nodal  displacements 
parameters  in  the  Newmark  integration  scheme 
strain  components  (1  =  1,2,. ..,6) 
orientation  of  m-th  layer  (m  =  1,2,...,L) 
density  of  m-th  layer  (m  =  1,2,...,L) 
stress  components  (1  *  1,2,..., 6) 

finite-element  Interpolation  functions  (1  =  1,2 . n) 

bending  slope  (rotation)  functions 
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Introduction 

The  transient  behavior  of  isotropic  and  homogeneous  plates  has  long 
been  a  subject  of  interest  {see  Lamb1  for  thin  plates  and  David  and 
Lawhead2  for  thick  plates).  For  many  years,  the  classical  (Polsson- 
Klrchoff)  plate  theory  (CPT),  in  which  normals  to  the  midsurface  before 
deformation  are  assumed  to  remain  straight  and  normal  to  the  midsurface 
after  deformation  (l.e.,  transverse  shear  strains  are  zero),  has  been 
used  to  calculate  frequencies,  static  response,  and  dynamic  response 
under  applied  loads.  Recent  studies  in  the  analysis  of  plates  have 
shown  that  the  effect  of  the  transverse  shear  strains  on  the  static  and 
dynamic  response  of  plates  Is  significant.  For  example,  the  natural 
frequencies  of  vibration  predicted  by  the  classical  plate  theory  are  25% 
higher,  for  plate  slde-to-thlckness  ratio  of  10,  than  those  predicted  by 
a  shear  deformation  theory  (SOT).  In  transient  analysis  of  plates  the 
classical  plate  theory  predicts  unrealistically  large  phase  velocities 
In  the  plate  for  shorter  wave  lengths.  The  Timoshenko  beam  theory3, 
which  includes  transverse  shear  and  rotary  Inertia  effects,  has  been 
extended  to  Isotropic  plates  by  Re1ssner^*5  and  MindllnS,  and  to 
laminated  anisotropic  plates  by  Yang,  Norris,  and  Stavsky7.  A 
general izatlon  of  the  von  (Carman  nonlinear  plate  theory  for  Isotropic 
plates  to  Include  the  effects  of  transverse  shear  and  rotary  Inertia  In 
the  theory  of  orthotropic  plates  Is  due  to  Medwadawskl®,  and  that  for 
anisotropic  plates  Is  due  to  Ebcloglu9. 

With  the  Increased  application  of  advanced  fiber  composite 
materials  to  jet  engine  fan  or  compressor  blades,  and  In  high 
performance  aircraft,  studies  Involving  transient  response  of  plates 
made  of  such  materials  are  needed  to  assess  the  ability  of  these 


materials  to  withstand  the  forces  of  impact  due  to  foreign  objects 
(e.g.,  the  Ingestion  of  stones,  nuts,  and  bolts,  hailstones,  or  birds  In 
jet  engines).  Previous  Investigations  Into  the  linear  transient 
analysis  of  composite  plates  Include  Moon's111,11  Investigation  of  the 
response  of  Infinite  laminated  plates  subjected  to  transverse  Impact 
loads  at  the  center  of  the  plate;  Chow's1 ^  study  of  laminated  plates 
(with  transverse  shear  and  rotary  inertia)  using  the  Laplace  transform 
technique;  Wang,  Chou,  and  Rose's1-1  investigation,  by  the  method  of 
characteristics,  of  unsymmetrical  orthotropic  laminated  plates;  and  Sun 
and  Whitney's1 4*1  ^  study  of  plates  under  cylindrical  bending.  More 
recently,  the  present  author1 ® >17  Investigated  the  linear  transient 
response  of  layered  anisotropic  composite  rectangular  plates  and 
presented  extensive  numerical  results  for  center  deflection  and 
stresses. 

Geometrically  nonlinear  transient  analysis  of  isotropic  plates  was 
considered  by  Hinton,  et  alJ8-20  an(j  ^^21  #  Hinton,  et  al.  used  the 
Mlndlln  element  while  Akay  used  a  mixed  finite  element  in  their  works. 

As  far  as  the  nonlinear  (geometric)  analysis  of  layered  anisotropic 
composite  plates  Is  concerned,  there  exist  no  previously  reported 
results  In  the  open  literature.  The  present  Investigation  Is  concerned 
with  the  geometrically  nonlinear  transient  analysis  of  layered  composite 
plates  under  applied  transverse  loads.  The  shear  deformable  element 
developed  by  the  author17  for  the  transient  analysis  of  layered 
composite  plates  Is  modified  to  Include  the  nonlinear  strain- 
displacement  relations  of  the  von  Karman  plate  theory.  Numerical 
results  are  presented  to  show  the  effect  of  lamination  scheme,  plate 
thickness,  nonlinearity,  boundary  conditions,  and  loading  on  the 
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transient  response  of  plates.  The  numerical  results  included  here  for 
layered  composite  plates  are  not  available  in  the  literature,  and 
therefore  should  be  of  interest  to  designers  of  composite-plate 
structures  and  numerical  analysts  and  experimentalists  In  evaluating 
their  techniques. 

Review  of  the  Equations  of  Motion  of  Anisotropic  Plates 
The  nonlinear  theory  of  laminated  anisotropic  plates  to  be  reviewed 
is  based  on  the  combination  of  the  Timoshenko-type  theory  and  the  von 
Karman  plate  theory.  The  theory  is  known  to  be  able  to  predict 
accurately  the  global  behavior.  However,  It  Is  not  accurate  enough  to 
predict  edge  stresses  and  hence  delamination.  The  theory  assumes  that 
the  stresses  normal  to  the  midplane  of  the  plate  are  negligible  when 
compared  to  the  inplane  stresses,  and  normals  to  the  plate  midsurface 
before  deformation  remain  straight  but  not  necessarily  normal  to  the 
midsurface  after  deformation. 

Equations  of  motion 

The  plate  under  consideration  is  composed  of  a  finite  number  of 
orthotropic  layers  of  uniform  thickness,  with  principal  axes  of 
elasticity  oriented  arbitrarily  with  respect  to  the  plate  axes.  The  x 
and  y-coordlnates  of  the  plate  are  taken  In  the  midplane  (q)  of  the 
plate.  As  in  all  Timoshenko-type  theories,  the  displacement  field  is 
assumed  to  be  of  the  form 

^(x.y.z.t)  =  u(x,y,t)  +  z<|,x(x,y,t) , 

u2(x,y,z,t)  -  v(x,y,t)  +  z  <|/y(x,y, t),  (1) 


u3(x,y,z,t)  »  w(x,y,t). 


■ 
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Here  t  is  the  time;  ult  113,  U3  are  the  displacements  in  the  x,  y,  z 

directions,  respectively;  u,  v,  w  are  the  associated  midplane 

displacements;  and  4.  and  4.  are  the  slopes  in  the  xz  and  yz  planes  due 

x  y 

to  bending  only.  The  strains  in  the  von  Karman  plate  theory  (which 
accounts  for  moderately  large  deflections  and  small  strains)  can  be 
expressed  in  the  form 


_  au  1  .aw.  2  A  ,  a' ^x_  _  0 
ei  ~  ax  +  2  (ax)  +  z  ax  =  ei  +  2,ci 


_  -  QL  x  A.  (0T-\C  4.  7 

e2  "  ay  2  'ay'  2  ay  _  fc2  T  ‘*2 

,  ,  au  a_v  aw  aw  , _  x.  V  =  0 

t-C  in  T  a  w  ▼  aw  aw  ▼  ^  \  a  w  T  T 


1  ,aw»2 


y  0 

=  Eo  +  ZlO 


ay  ax  ax  ay 


.  ,  aw 

*x  +  ax 


'ay 

c4  =  *y  + 


ax  •  -  ^6 

aw 


ay 


(?) 


wherein  the  squares  of  the  first  spatial  derivatives  of  u,  v, 

4;  and  4,  are  neglected.  The  strain  does  not  enter  the  equations 
x  y  j 

because  the  constitutive  relations,  to  be  given  shortly,  are  based  on 
the  plane-stress  assumption.  Note  that  the  transverse  shear 
strains,  and  e,.,  are  constant  through  the  thickness.  If  a  linear 
distribution  of  the  transverse  shear  strains  through  the  thickness  is 
desired,  one  must  add  higher  order  terms  in  z  to  the  displacements  uj 
and  U2,  and/or  U3:  with  each  additional  term,  an  additional  dependent 
variable  Is  introduced  Into  the  problem. 

Neglecting  the  body  moments  and  surface  shearing  forces,  we  write 
the  equations  of  motion  in  the  presence  of  applied  transverse  forces,  q. 


as 
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A 
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h 
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"l.x  +  N6,y  =  Pu,tt  +  R*x,tt 


N6,x  +  N2,y  ”  Pv,tt  +  R<ly,tt 

°l,x  +  °2,y  =  +  q(x’y>t)  +  W("*V 


(3) 


M,  +  M,  -  0  =  Id-  +  Ru 
l.x  6,y  1  vx,tt  ,tt 


M_  +  -  0,  =  let-  **  +  Rv 

6,x  2,y  2  y,tt  ,tt 


where  P,  R  and  I  are  the  normal,  coupled  normal -rotary,  and  rotary 
inertia  coefficients. 


^/2  n  m+1  n  f  , 

(P.R.r)  =  /  (l,z,z^)pdz  =  l  /  (l,z,zZ)p'm'dz  (4) 

-h/2  m  z 

m 

being  the  material  density  of  the  m-th  layer,  N-j ,  Qi ,  and  M-|  are 
the  stress  and  moment  resultants  defined  by 


h/2  h/2 

,M1)  =  /  (l,z)a1dz  ,  (QpQo)^  (o5,a4)dz,  (5) 

-h/2  c  -h/2 


and  A(w,Ni)  the  contribution  due  to  the  nonlinear  terms. 


aw  BNi  5Nfi  aw 

■  f?  <sr  *  »r)  ♦  i  < 


aN.  aN- 

_ 2.  .  _ £\ 

ay  vax  ay  ' 


(6) 


Here  (i  =  1,2,. ..,6)  denote  the  in-plane  stress  components 

(°1  *  V  ■  V  °4  •  v  °s  *  °xz  and  °6  ■  °xy>- 

If  one  plane  of  elastic  symmetry  parallel  to  the  plane  of  each 

layer  exists,  the  constitutive  equations  for  the  plate  can  be  written  In 

the  form  (see  Jones22) , 


N1 

a 

- 1 

C-*. 

CD 

_ 1 

0 

Ej 

l  °2 

"'a* 

e4 

Mi 

-Bj1  D1j  - 

*  U 

cn 

A55  - 

e5 

The  Atj,  Bfj,  Dij  (1,j  =  1,2,6),  and  *ij  (i ,j  =  4,5)  are  the  inplane, 
bendlng-lnplane  coupling,  bending  or  twisting,  and  thickness-shear 
stiffnesses,  respectively: 


Zm+1  0  zm+l 


u  -l  /.  k1kj°i j >dz 

(8) 


m  z. 


m 


m  z. 


Here  z„,  denotes  the  distance  from  the  mid-plane  to  the  lower  surface  of 
the  m-th  layer,  and  k-|  are  the  shear  correction  coefficients. 

Variational  Formulation 


The  variational  form  of  Eqs.  (3)  and  (7)  is  given  by 


/  {6u(Pu  + 

Q 

*x, 

,tt^ 

*  6u,xNl  < 

■  6U,A  *  6,(p,,u 

+  %,«> 

+  6v  xNg  +  6v 

,yN2 

+  6w(Pw  tt)  + 

6w  Q,  +  6w 

»A  1  5 

,y°2 

.  asw  ,aw 
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asw  ,aw 
ay  'ax 

N6  *  1?  N2> 

+  ^x^x.tt 

+  Ru 

,tt^ 

+  6cp  M 
x,x  1 

+  6*x,A  * 

*xQl 

+  V^y.tt 

+  Rv 

,tt^ 

+  6(l» 

y.x  6 

*  ‘VA  * 

64>yQ2 

+  6wq}dxdy 

+  /  (6unNn  +  6usNns)dS  +  /  6wVds  +  /  (6^  +  64-^^ds, 

(9) 


wherein  Ni ,  and  are  given  In  terms  of  the  generalized 
displacements  by  Eq.  (7),  and  V,  N„  and  Nns,  and  and  MnS  are  the 
shear  force,  normal  and  tangential  in-plane  forces,  and  normal  and 
twisting  bending  moments  defined  on  the  boundary  r,  respectively: 
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N  =  n  n  N,  +  ?n  n  N£  +  n  n  N, 
n  xxl  x  y  6  y  y  2 

N  =  n  n  (N„  -  N. )  +  (n  n  -  n  n  )N, 
ns  x  y'  2  1  x  x  y  y  6 

M  =  n  n  N,  +  2n  (i  N,  t  n  i)  H 
n  xxl  x  y  6  y  y  2 

Mns  ■  nxny(H2  ‘  Ml>  +  (nx\  -  Vy,M6’  (ln) 

A 

where  n  =  (nx,n y)  is  the  unit  vector  normal  to  the  boundary  r.  The 
variational  formulation  indicates  that  the  essential  (i.e.,  geometric) 
and  natural  boundary  conditions  of  the  problem  are  given  by: 

essential  specify:  un,  us,  w,  <j,n,  <|,s 

natural  specify:  Nn,  %s.  V,  M„,  fV»s*  (H) 

wherein  un  and  us,  for  example,  denote  the  normal  and  tangential 
components  of  the  Inplane  displacement  vector,  u  =  (u,v). 

A  Shear-Deformable  Finite-Element  for  Laminated  Plates 
Here  we  present  a  finite-element  model  ;<$soc1ated  with  the 
nonlinear  equations  governing  the  motion  of  layered  composite  plates. 

The  element  is  an  extension  of  the  penalty  plate-bending  element 
developed  for  the  linear  static  and  dynamic  analysis  of  layered 
composite  plates  by  the  present  author-! 7,23. 

Consider  a  finite-element  analog,  of  the  midplane  of  the  plate, 
Q.  Over  a  typical  element,  of  the  mesh  Qh,  each  generalized 
displacement  U  Is  Interpolated  spatially  by  an  expression  of  the  form, 

U  *  I  lUt)*.(x,y) 

1  1  1 


(12) 
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where  U-j  is  the  value  of  U  at  node  i  at  time  t,  4^  is  the  finite-element 

interpolation  function  at  node  1  and  r  is  the  number  of  nodes  in  the 

element.  For  the  sake  of  simplicity,  we  use  the  same  interpolation 

function  for  each  of  the  generalized  displacements,  (u,v,w,<p  ,4.  ). 

x  y 

Specifically,  the  nine-node  isoparametric  rectangular  element  with  five 
degrees  of  freedom  is  used  In  the  present  study. 

Substituting  Eqs.  (7)  and  (12)  into  Eq.  (9),  we  obtain  the 
following  equation: 

[K]{a}  +  =  {F}  (13) 

Here  {a}  is  the  column  vector  of  the  nodal  values  of  the  generalized 
displacements,  [K]  is  the  matrix  of  stiffness  coefficients,  [M]  is  the 
matrix  of  mass  coefficients,  and  {F}  is  the  column  vector  containing  the 
boundary  and  transverse  force  contributions.  The  elements  of  [K]  and 
[M]  are  given  In  Appendix  A. 

It  should  be  pointed  out  that  the  element  stiffness  matrix  [K]  is 
nonlinear  and  unsymmetrlc  in  the  present  formulation  (see  Appendix  A), 
the  nonlinearity  being  due  to  the  nonlinear  terms  appearing  in  the 
variational  formulation.  The  nonsymmetric  nature  of  the  stiffness 
matrix  needs  some  explanation.  To  this  end  consider  the  following  terms 
from  the  variational  formulation, 


+  A 


36u  rau 
11  ax  '-ax 


1  .aw.2_ 

+  ?  (57)  ]  + 


asw  aw  _au  1  .aw.  2- 
11  ax  ax  Lax  +  7  ^ax'  -* 


The  first  term  contributes  to  the  stiffness  coefficients  and 
(see  Appendix  A),  whereas  the  second  term  contributes  to  the  stiffness 
coefficients  K^j  and  K^j.  Note  that 


43 


■  / 


41  k7 


aw. 


a<i» 


ax  ax' 


J  dxdy  + 


are  not  the  same. 


To  complete  the  discretization,  we  must  now  approximate  the  time 
derivatives  appearing  In  Eq.  (13)  Here  we  use  the  Newmark  direct 
integration  method^,  with  a  =  0.5  and  p  =  0.25  (corresponding  to  the 
constant-average-acceleration  method).  The  scheme,  although 
unconditionally  stable  for  linear  problems,  is  not  proved  to  be  stable 
for  nonlinear  problems. 

Equation  (13)  can  be  expressed,  after  the  application  of  Newmark 's 
integration  scheme,  in  the  form 

-  <Fln>nt,  »*) 

where 

[K]  =  [K]  +  a0[M],  {F}  =  {F}n+1  +  [M](ao{A}n  +  a^A},  +  a2{A}n). 

=  1  /  ( p At ^ )  ,  flj  =  a^At  »  a2  =  ^  •  (15) 

Once  the  solution  is  known  at  tn+1  =  (n+l)&t,  the  first  and  second 
derivatives  (velocity  and  accelerations)  of  {a}  at  tn+]  can  be  computed 
from 

Wn*l  ■  ao«4tn+l  -  WB1  -  al<‘>n  -  a2Wn 

“  Wn  *  a3,4tn  +  (16) 

where  83  =  (1  -  a)At,  and  84  =  aAt. 


J 
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All  of  the  operations  Indicated  above  (except  for  eqn.  (16))  can  be 
performed  at  the  element  level,  and  the  assembled  form  of  eqn.  (14)  can 
be  obtained  for  the  whole  problem: 

[K]  =  [KL  +  KN({Mn+1)]  (17) 

where  [KL1  denotes  the  linear  stiffness  matrix  and  denotes  the 
nonlinear  (geometric)  stiffness  matrix.  Because  [KN]  depends  on  the 
unknown  solution  {A}n+^  ,  the  assembled  equation  must  be  solved 
iteratively  until  a  convergence  criterion  is  satisfied  at  time,  t  = 
tn+i.  In  the  present  study  the  Picard  type  successive  iteration  scheme 
is  employed.  In  this  scheme,  the  nonlinear  equations  (14)  are 
approximated  by  the  following  equation 

(18) 

where  r  is  the  iteration  number.  In  other  words,  the  nonlinear 
stiffness  matrix  for  the  (r  +  1 )-th  iteration  is  computed  using  the 
solution  vector  from  the  r-th  iteration.  Such  successive  iterations  are 
continued  until  the  error 

E  =  [  z  |  aJ  -  aJ+1  I2  e|a!J|2]  (19) 

i  =1  1  1  1  1 

for  any  fixed  time  t  *  tn+j ,  is  less  than  or  equal  to  some  preassigned 
value  (say,  one  percent  or  less).  The  iteration  at  time  t  =  tn+l  is 
started  by  using  the  converged  solution  at  t  >  tn  (at  t  =  to  =  0,  the 
initial  conditions  are  assumed  to  be  known). 


Numerical  Results  and  Discussion 


In  the  present  study  the  nine-node  rectangular  isoparametric 

element  was  employed.  The  element  has  either  three  (w,<px  ,(|>y)  or 

five  (u,v,w,<j>  ,4>  )  degrees  of  freedom  (DOF)  per  node.  Because  the 
x  y 

element  accounts  for  the  transverse  shear  strains,  reduced  integration 
was  employed  to  evaluate  the  shear  terms  numerically.  In  other  words, 
the  2x2  Gauss  rule  was  used  to  Integrate  the  shear  energy  terms  and  the 
3x3  Gauss  rule  was  used  to  integrate  the  bending  and  Inertia  terms.  All 
computations  were  done  In  double  precision  on  an  IBM  3032.  Due  to  the 
biaxial  symmetry  of  the  problems  discussed,  only  one  quadrant  of  the 
plate  was  analyzed.  In  all  of  the  ninerical  examples  considered  here, 
zero  Initial  conditions  are  assumed  and  damping  is  neglected. 

Because  no  estimate  on  the  time  step  for  the  nonlinear  analysis  is 
available  in  the  literature,  the  critical  time  step  of  a  conditionally 
stable  finite  difference  scheme  was  used  as  a  starting  time  step,  and  a 
convergence  study  was  conducted  to  select  a  time  step  that  yielded  a 
stable  and  accurate  solution  while  keeping  the  computational  time  to  a 
minimum.  The  following  two  estimates  were  used  in  the  present  study: 

Atj  <  0.25  (phD)1/2(Ax)2  (20) 

2 

At2  <  [p(l  -  v2) E/ { 2  +  11  -  v)  ^  [1  +  1.5  (^)2}]1/2AX  (21) 

Here  D  =  Eh3/[i2(l  -  v2)],  and  ax  is  the  minimum  distance  between  the 
element  node  points.  Estimate  (20),  due  to  Leech^5,  was  derived  for 
thin  plates,  and  estimate  (21),  due  to  Tsui  and  Pin  Tong26,  was  derived 
for  thick  plates. 

First,  in  order  to  prove  the  validity  of  the  present  formulation, 
the  results  obtained  in  the  present  study  for  isotropic  plates  are 
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compared  with  those  available  in  the  literature: 

1.  Simply  supported  rectangular  plate  under  suddenly  applied  patch 
loading  (see  Fig.  1  and  Ref.  [27]) 
a  =  /?,  b  =  1 ,  h  =  0.1  and  0.2,  At  =  0.1 ,  q  =  qQH(t) 


,10“2,  0  <  x,y  <  0.2 
0  ,  x,y  >  0.2 


(22) 


E^  =  Ep  =  1.0,  Vjp  —  0.3,  p  =  1.0. 

2.  Simply  supported  square  plates  under  suddenly  applied  uniformly 
distributed  loading  (see  Fig.  2  and  Ref.  [21]) 

-4  2 

a=b  =  243.8  cm,  h  =  0.635  cm.  At  =  0.005  sec.,  qQ  =  4.882  x  10  N/cm  , 
Ej=E2  =  7.031 xl 05  N/cm2,  v12  =  0.25,  p  =  2.547  x  10"6  N  sec2/cm4. 

(23) 

3.  Clamped  circular  plate  under  suddenly  applied  uniformly 
distributed  loading  (see  Fig.  3  and  Ref.  [19]) 

R(radlus)  =  100  In.,  h  =  1 0  and  20  In.,  At  =  2.5  sec.,  qQ  =  1.0  psi , 

Ei  =  E2  =  100  psi,  v12  =  0.3,  p  =  10  lb  sec.2/1nch4.  (24) 

1.  Simply  supported,  rectangular,  isotropic  and  orthotropic 
plates.  A  4x4  (nonuniform)  mesh  Is  used  in  the  quarter  plate  with  3  DOF 
per  node.  The  center  deflection  and  bending  moments  of  the  present 
linear  analysis  are  compared  with  the  analytical  thick-plate  and  thin- 
plate  solutions  of  Reismann  and  Lee27  in  Fig.  4.  we  note  pronounced 
difference  between  the  solutions  of  the  two  theories.  The  deviation 
between  the  solutions  of  the  shear  deformation  theory  and  the  classical 
theory  increases  with  plate  thickness.  The  present  finite  element 
solutions  for  the  center  deflection  and  bending  moment  are  in  excellent 
agreement  with  those  of  Reismann  and  Lee^. 


Figure  1.  Finite-element  mesh  and  boundary  conditions 

for  isotropic  rectangular  plates  under  suddenly 
applied  pressure  loading  at  the  center. 


v-w-V0 


BC1:  u(0,y)«v(x,0)«0 
BC2:  u(y,0)"v(0,y)«0 


Figure  2.  Finite-element  mesh  and  boundary  conditions 
for  isotropic  square  plates  under  suddenly 
applied  pressure  loading. 


r 


Nonil lmensional ized  center  bending  moment,  Nond lmensional lzed  deflection 


The  center  deflections  obtained  In  the  present  linear  and  nonlinear 
analysis  are  presented  In  Fig.  5.  To  show  the  effect  of  material 
orthotropy  on  the  center  deflection,  results  of  the  linear  analysis  of 
an  orthotropic  plate  (Ej/Eg  =  25,  v12  =  0.25)  are  also  included  in  Fig. 

5.  The  effect  of  orthotropy  and  geometric  nonlinearity  on  the  amplitude 
and  period  of  nondlmenslonal 1  zed  center  deflection  of  orthotropic  plates 
(El/E2  s  25,  v12  =  0.25)  Is  apparent  from  the  plots  presented  In  Fig.  5; 
the  effect  of  orthotropy  is  to  decrease  both  amplitude  and  period  of  the 
center  deflection,  and  effect  of  the  nonlinearity  is  to  decrease  the 
amplitude,  and  smoothen  the  solution  somewhat.  The  effect  of  the  plate 
thickness  on  the  amplitude  and  period  of  the  nondimensionalized  center 
deflection  (linear)  w  =  10wE2h3/qQh4,  can  be  seen  from  the  plots  in  Fig. 

6.  The  amplitude  decreased  about  30%  for  a  decrease  in  thickness 
from  £  ■  0.2  to  ^  =  0.1 . 

2.  Simply-supported  square  plates  (BC1).  A  2x2  (uniform)  mesh  Is 
used  In  the  quarter  plate  with  5  DOF  per  node.  Results  of  the  present 
nonlinear  analysis  (see  Fig.  7)  agree  closely  with  the  mixed  finite 
element  results  of  Akay21.  The  plots  of  center  deflection  versus  time 
for  various  loads  are  shown  In  Fig.  7  along  with  the  load-deflection 
curve.  Note  that  amplitude  and  period  of  the  center  deflection  decrase 
with  Increasing  values  of  the  load.  Further  note  that  the  negative  peak 
also  Increases  with  decreasing  load. 

3.  _ Clamped  circular  plates.  Results  of  the  linear  transient 

anlysls  of  the  clamped  square  plate,  for  two  different  thicknesses,  are 
presented  In  Figs.  8a  and  b.  The  present  results  agree  with  those 
reported  by  H1ntonl9.  The  effect  of  the  decrease  in  thickness  is  to 
Increase  the  amplitude  and  period  of  the  center  deflection  and  stress. 


Center  deflection  (wxl(r) 
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Figure  5. 


Center  deflection  versus  time  for  isotropic  (see  eqn.  (22)) 
and  orthotropic  (E,=25,  E2=l  >  v,2=0.25,  G, 2=(H3=®23=0’ 5E2^ 
rectangular  plates  under  suddenly  applied  Hodfng  at  the 


center  square  area(see  Figure  1). 


Nondimensionalized  center  deflection 


20 


Time,  t 


Figure  6.  Effect  of  plate  thickness  on  the  transient  response  of 

orthotropic  plates  (see  Figure  5  for  material  properties) 
under  suddenly  applied  patch  loading  (i.e.,  load  on  the 
square  area  at  the  center  of  the  rectangular  plate). 


Center  deflection,  w  (cm) 
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Figure  7.  Nonlinear  transient  response  of  isotropic  plates 

(see  eqn.  (23)  and  Figure  2  for  the  data  and  finite 
element  mesh  and  boundary  conditions)under  suddenly 
applied  uniformly  distributed  loads. 
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(a)  Clamped  circular  plate,  h  *  20". 


i  1 
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;  i 


(b)  Clamped  circular  plate,  h  *  10". 

Figure  8.  Transient  response  of  a  clamped,  isotropic,  circular 
plate  under  suddenly  applied  transverse  load(see 
Figure  3  for  the  finite-element  mesh  and  boundary 
conditions). 
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Figure  8a  also  contains  plots  of  center  deflection  of  the  linear  and 
nonlinear  analysis.  The  effect  of  the  (von  Karman  type)  nonlinearity  is 
to  decrease  the  amplitude  and  period  of  the  center  deflection  and 
stress. 

The  results  of  the  nonlinear  analysis  of  layered  composite  plates 
are  presented  next.  Figures  9  and  10  contain  plots  of  the  center 
deflection  versus  time  for  simply  supported  (BC1)  two-layer  cross-ply 
[0°/90°]  and  angle-ply  [45°/-45°]  plates,  respectively.  The  same  data 
(except  for  the  material  properties)  as  In  Eq.  (22)  is  employed.  The 
effect  of  coupling  between  the  inplane  displacements  and  bending 
deflections  is  to  increase  the  amplitude  of  the  center  deflection  in  the 
linear  analysis.  Figure  10  also  contains  the  results  of  uniformly 
loaded  angle-ply  plates.  It  Is  apparent  from  this  figure  that  as  the 
load  is  changed  from  uniformly  distributed  to  patch  loading,  the 
response  curves  change  from  a  regular  to  Irregular  curve  (and  for  point 
loading  the  response  becomes  spiky). 

The  plots  of  the  center  deflection  versus  time  for  single-layer 
(0°),  cross-ply  [0°/90o]f  and  angle-ply  [45°/-45°]  square  plates  under 
suddenly  applied  uniform  loading  are  presented  in  Figs.  11  and  12, 
respectively;  2x2  (uniform)  mesh  with  5  DOF  is  used.  The  boundary 
conditions  used  in  this  set  of  problems  (BC2)  are  different  from  the 
simply  supported  boundary  conditions  (BC1)  used  earlier;  the  two  types 
of  boundary  conditions  are  shown  in  Fig.  2.  The  data  used  Is  the  same 
as  that  in  Eq.  (23)  (except  for  the  material  properties).  The  effect  of 
the  nonlinearity  and  the  lamination  scheme  on  the  amplitude  and  period 
of  the  center  deflection  Is  apparent  from  the  results.  Note  that  the 
rate  of  Increase  In  the  amplitude  of  the  center  deflection  decreases 
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Figure  9.  Transient  response  of  two-layer  cross-ply  (0°/90°) 
rectangular  plates  under  suddenly  applied  patch 
loading  (see  eqn.(22)  and  Figure  1  for  the  data  and 
mesh  information;  see  Figure  5  for  the  material 
properties) . 


Nonlinear  transient  response  of  two-layer  angle-ply 
[45°/-45<-3  rectangular  plates  under  suddenly  applied 
patch  loading  (see  eqn.  (22)  and  Figure  1  for  data). 


Figure  11.  Nonlinear  transient  response  of  single-layer 
(0°)  orthotropic  square  plates  under  suddenly 
applied  pressure  loading  (see  eqn.  (23)  and 
Figure  2  for  the  data  and  mesh);  BC2  boundary 


conditions  were  used. 
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(b)  Angle-ply  [45°/-45°]  plates 

Figure  12.  Nonlinear  transient  response  of  laminated  square  plates 
under  suddenly  applied  uniform  loading  (see  Figure  10 
for  other  information). 
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with  increasing  load.  Further  note  that  the  maximum  center  deflection 
for  the  cross-ply  plates  is  about  10%  larger  than  that  of  the  angle-ply 
[45°/-45°]  plates. 

Lastly,  the  effect  of  boundary  conditions  (BC1  and  BC2)  on  the 
center  deflection  of  two-layer  square  plates  under  uniformly  distributed 
load  was  investigated;  Figure  13  contains  results  of  the 
investigation.  It  is  clear  that  the  response  of  cross-ply  as  well  as 
angle-ply  plates  under  BC1  differ  significantly  from  the  results  of 
plates  under  BC2.  Note  also  that  suppression  of  the  inplane  degrees 
freedom  results  In  about  50%  decrease  in  the  center  deflection  of  angle- 
ply  [45°/-45°]  plates. 


Summary  and  Conclusions 

A  shear  flexible  finite  element  that  accounts  for  the  von  Karman 
strains  is  employed  in  the  transient  analysis  of  layered  composite 
plates.  The  present  results  for  isotropic  plates  are  very  close  to  the 
analytical  and  other  finite  element  solutions  available  in  the 
literature.  The  present  results  of  layered  composite  plates  should 
serve  as  reference  solutions  for  future  investigations.  Although  the 
paper  contains  results  for  certain  geometries,  loadings,  lamination 
scheme,  and  material  properties,  it  should  be  pointed  out  that  the 
element  developed  herein  can  be  employed  to  laminated  plates  of 
arbitrary  geometry,  lamination  scheme,  material  properties,  boundary 
conditions,  and  loading  (only  limitations  are  those  implied  In  the 
laminated  plate  theory  used).  The  present  analysis  does  not  account  for 
material  damping  effects;  studies  of  composite  plates  accounting  for 
material  damping  and  material  nonlinearities  are  awaiting  attention. 


Center  deflection,  w(cm) 


Time  (txlO2) 


Figure  13.  Effect  of  the  boundary  conditions  (BC1  and  BC2)  on  the 
nonlinear  t  nsient  response  of  laminated  square  plates 
under  suddenly  applied  transverse,  uniformly  distributed 
loading  (see  eqn.  (23)  for  the  data  and  Figure  2  for  the 
mesh  and  boundary  conditions;  orthotropic  properties  of 
layers  are  given  in  Figure  11). 
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APPENDIX  A  :  ELEMENTS  OF  STIFFNESS  AND  MASS  MATRICES 
stiffness  matrix: 


[KU] 

CK12] 

CK13] 

[k143 

CK15! 

CK21] 

[K22] 

CK23! 

CK24! 

CK25! 

[K]  = 

[k311 

C«32.l 

CK33] 

CK?4! 

CK?5! 

+  [k|31 

+  CK34] 

+  Ck|5! 

[K41] 

CK42] 

[<l3] 

CK44! 

CK45] 

+  CK43] 

[K51] 

[kb2] 

Ck”] 

CK54! 

CK55! 

+  [<23] 

mass  matrix 


PCS! 

[0] 

[0] 

RCS] 

[0] 

Cn] 

PCS! 

[0] 

[0] 

RCS] 

[0] 

CO] 

PCS] 

CO] 

CO] 

RCS! 

CO! 

Co] 

i[s] 

[0] 

[0] 

RCS! 

[0] 

[0] 

ICS] 

The  matrix  coefficients  are  given  by 

CK11]  =  An[SXX]  +  Alfi([SXy]  +  [SXy]T)  +  A66[Syy], 

[K12]  =  A12[SXy]  +  A16[SXX]  +  A^S™]  +  A66[SXy]T 

-  C<211T, 

[K13]  =  An[RXX]  +  A12[RXy]  +  Alfi([RXy]  +  [RXy]T 

♦  [Ryx])  +  AggC^]  +  A66(CRyy]T  +  [Ryy])  =i[K313T, 


33 


[K14]  =  B11[SXX]  +  B16(CSXy]  +  [Sxy]T)  +  B66[Syy]  =  [K41]T, 

[K15]  =  B12[SXy]  +  B16[SXX]  +  B2g[Syy]  +  B66[SXy]T 

[KZZ]  -  A22[Syy]  +  A26([SXy]T  +  [SXy])  +  A66[SXX], 

[K23]  =  A12[Rxy]T  +  A22[Rjy]  +  A26([RxyHT  +  [R^y] 
+l[Rf])+A16rRXX1  +  A66([RX/]T+[R;X]) 

=  \  [k321t  , 

[K24]  =  B12fSXy1T  +  B26CSyyl  +  B16CSXX7  +  866rsxy_7 

-  [K421T, 

[K25]  =  B^c^]  +  B26([SXy]  +  CSXy]T)  +  B66[SXX] 

«  [K52]T, 

[Kf]  -  A55[SXX]  +  A45([SXy]  ♦  [SXy]T)  ♦  A44[Syy], 

.  _  8$.  8$.  _  8<i>.  8<t>.  9<t> .  3<j> . 

£K2  3  =  7  /Re^Nl  ax"  ax1  +  N6  W~  Ty~  +  VyL"5xJ  } 
a<i).  a*. 

+  ^2  W  ^T]dxdy’ 

[K?4]  -  A55rsx0]  +  A45[Sy0]  =  [Kf  1T, 

[K34]  =  Bu[Rjxl  +  B16([Rxy]  +  [RXy  ]T  +  [RXX]) 

+  B12[RXy]T  +  B26CRyy7  ♦  B66(f<y7  ♦  [RXyD 


+ 


o<p  -  A45rsxoi  t  A44rsyo] .  pep'. 

[Kf]  -  B^J  ♦  B22[Rjy]  *  B26([Rjy]T 

CRf ]  +  [Ryy])  +  B66([R*y]T  +  [Rxx])  =  2[K53]T 
[K44]  =  DU[SXX]  +  ^^CSXy3  +  [SXy]T)  +  D6g[Syy]  + 

+  A55[S3‘ 

[K451  -  D12[SXy3  +  D16[SXX]  +  D26[Syy]  +  D66[SXy]T  + 

+  a45[s]  =  [k54]t, 

[K55]  =  D26([SXyl  +  [SXy]T)  +  Dgg[SXX]  +  D22[Syy]  + 


S1j  =  /  W  dr^  dxfiy  ’  *-,r]  =  S?j  5  Sij» 


Rc  1  7  ^a^  bt,'  aiT  dxdy  ’  =  °»x»y» 


Jf  -  2  +  a  f^w*2  +  ?A  5w 

P  A11  W  A12  ("ay)  2Aifi  -^7  ♦ 


16  3x  ay 


^2  =  A12  (S)2  +  A22  +  2A 


3w,2  3w  aw 

26  ax  ay  » 


m  _*  r£iL\  2  ,  »  /.dW\  2  ,  -«  aw  3w 

"6  *16  'ax>  +  a26  ^ay)  2ARK  to1  *7  • 


66  ax  ay 
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